---
title: "Fixed Imaging Quantification"
output:
  html_document:
   code_folding: hide
---

```{r toolbox/importation, include=FALSE, warning = FALSE}
library(ggplot2)
library(moments)
library(dplyr)
library(emmeans)
library(ggsignif)
library(MASS)
library(lme4)
library(lmtest)
library(RColorBrewer)

ens0 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/older models/Fixed ENS Image Quant Combined Data v5.csv", header=TRUE, stringsAsFactors=FALSE)

cort0 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/older models/Fixed cortex Image Quant Combined Data v2.csv", header = TRUE, stringsAsFactors = FALSE)

ens_tubb0 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/B3T Avg Intensity/ens tubb3 v3.csv", header=TRUE, stringsAsFactors=FALSE)

cort_tubb0 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/B3T Avg Intensity/cortex tubb3 v1.csv", header = TRUE, stringsAsFactors = FALSE)

ens0$Groups.f <- factor(ens0$group, levels = c("Control","But", "LPS", "But+LPS",
                                             "PFF", "PFF+But", "PFF+LPS", "PFF+But+LPS",
                                             "2XPFF", "4XPFF", "4XPFF+But", "6XPFF"))
ens0$attempt.no <- factor(ens0$attempt)

ens0$pff.qual <- ifelse(ens0$attempt.no == 6 & 
                            (ens0$group == "6XPFF" | 
                             ens0$group == "6XPFF+LPS" | 
                             ens0$group == "6XPFF+But" | 
                             ens0$group == "6XPFF+But+LPS"), 
                       1,0)



ens0$but.qual <- ifelse(ens0$attempt.no == 19 & (ens0$Groups.f == "But" | 
                                                  ens0$Groups.f == "But+LPS" | 
                                                  ens0$Groups.f == "PFF+But" | 
                                                  ens0$Groups.f == "4XPFF+But"), 1,0)

ens0$w.id <- factor(paste(ens0$attempt.no,"-",ens0$Groups.f,"-",ens0$well))

cort0$Groups.f <- factor(cort0$group, levels = c("Control","PFF","4XPFF"))
cort0$attempt.no <- factor(cort0$attempt)

ens_tubb0$attempt.no <- factor(ens_tubb0$attempt)

ens_tubb0$Groups.f <- factor(ens_tubb0$group, levels = c("Control","But", "LPS", "But+LPS",
                                             "PFF", "PFF+But", "PFF+LPS", "PFF+But+LPS",
                                             "2XPFF", "4XPFF", "4XPFF+But", "6XPFF"))

ens_tubb0$pff.qual <- ifelse(ens_tubb0$attempt.no == 6 & 
                            (ens_tubb0$group == "6XPFF" | 
                             ens_tubb0$group == "6XPFF+LPS" | 
                             ens_tubb0$group == "6XPFF+But" | 
                             ens_tubb0$group == "6XPFF+But+LPS"), 
                       1,0)

cort_tubb0$Groups.f <- factor(cort_tubb0$group, levels = c("Control","PFF","4XPFF"))
cort_tubb0$attempt.no <- factor(cort_tubb0$rat)
```

```{r unify dapi count cortex,  warning = FALSE}
cortdap.10 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/DAPI Counting Analysis/Cortex Quant (except 14)/MyExpt_Nuclei.edit.csv", header = TRUE, stringsAsFactors = FALSE)

cortdap.20 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/DAPI Counting Analysis/Cortex (Att 14)/MyExpt_Nuclei.csv", header = TRUE, stringsAsFactors = FALSE)

cortdap.11 <- cbind(cortdap.10[,1],cortdap.10[,2],cortdap.10[,7],cortdap.10[,10],cortdap.10[,12],cortdap.10[,9])
cortdap.21 <- cbind(cortdap.20[,1],cortdap.20[,2],cortdap.20[,8],cortdap.20[,11],cortdap.20[,13],cortdap.20[,10])

cortdap <- rbind(cortdap.11,cortdap.21)
cortdap <- data.frame(cortdap)

colnames(cortdap) <- c("image.no",
                       "nuclei.no",
                       "group",
                       "rat.no",
                       "well",
                       "im.no")

cortdap$id.maker <- paste(cortdap$rat.no,"-",
                          cortdap$group,"-",
                          cortdap$well,"-",
                          cortdap$im.no)

cortdap.short <- data.frame(aggregate(data = cortdap,
                          nuclei.no ~ id.maker,
                          function(nuclei.no) length(unique(nuclei.no))))

cort0$id.maker <- paste(cort0$attempt,"-",
                          cort0$group,"-",
                          cort0$well,"-",
                          cort0$image.no)

cort.m2 <- merge(cort0,cortdap.short, all = TRUE)
cort <- merge(cort0,cortdap.short,by="id.maker")
```

```{r unify dapi count ens,  warning = FALSE}
ensdap.10 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/DAPI Counting Analysis/ENS Quantification (except 14)/MyExpt_Nuclei.edit.csv", header = TRUE, stringsAsFactors = FALSE)

ensdap.20 <- read.csv("/Users/adambindas/Dropbox/Parkinson's disease ENS Paper/Data (Raw and Organized)/Immunostaining Image Quantification/DAPI Counting Analysis/ENS (Att 14)/MyExpt_Nuclei.csv", header = TRUE, stringsAsFactors = FALSE)

ensdap.11 <- cbind(ensdap.10[,1],ensdap.10[,2],ensdap.10[,7],ensdap.10[,10],ensdap.10[,12],ensdap.10[,9])
ensdap.21 <- cbind(ensdap.20[,1],ensdap.20[,2],ensdap.20[,8],ensdap.20[,11],ensdap.20[,13],ensdap.20[,10])

ensdap <- data.frame(rbind(ensdap.11,ensdap.21))

colnames(ensdap) <- c("image.no",
                       "nuclei.no",
                       "group",
                       "rat.no",
                       "well",
                       "im.no")

ensdap$id.maker <- paste(ensdap$rat.no,"-",
                          ensdap$group,"-",
                          ensdap$well,"-",
                          ensdap$im.no)

ensdap.short <- data.frame(aggregate(data = ensdap,
                          nuclei.no ~ id.maker,
                          function(nuclei.no) length(unique(nuclei.no))))

ens0$id.maker <- paste(ens0$attempt,"-",
                          ens0$group,"-",
                          ens0$well,"-",
                          ens0$image.no)

ens_tubb0$id.maker <- paste(ens0$attempt,"-",
                          ens0$group,"-",
                          ens0$well,"-",
                          ens0$image.no)

ens <- merge(ens0,ensdap.short,by="id.maker")
ens_tubb.merge <- merge(ens_tubb0,ensdap.short,by="id.maker")
```

```{r Initial Inspection/Organization,  warning = FALSE}
ens2 <- subset(ens0, 
                 quality.removal == 0 & sat.qual.test == 0 & 
                 pff.qual == 0 & but.qual == 0 &
                 attempt != 4 & attempt != 14 & attempt != 6)

ens3 <- subset(ens2,
                 X6xpff == 0 & X6xpff.but == 0 & X6xpff.lps == 0 & X6xpff.but.lps == 0 & 
                 X2xpff == 0 & X4xpff == 0 & X4xpff.but == 0)

ens4 <- subset(ens2,
                 X6xpff.but == 0 & X6xpff.lps == 0 & X6xpff.but.lps == 0)

ens2.nuc <- subset(ens, 
                 quality.removal == 0 & sat.qual.test == 0 & 
                 pff.qual == 0 & but.qual == 0 &
                 attempt != 4 & attempt != 14 & attempt != 6)

ens3.nuc <- subset(ens2.nuc,
                 X6xpff == 0 & X6xpff.but == 0 & X6xpff.lps == 0 & X6xpff.but.lps == 0 & 
                 X2xpff == 0 & X4xpff == 0 & X4xpff.but == 0)

ens4.nuc <- subset(ens2.nuc,
                 X6xpff.but == 0 & X6xpff.lps == 0 & X6xpff.but.lps == 0)


cort2 <- subset(cort0, quality.removal == 0 & sample.size == 0)
cort2.nuc <- subset(cort, quality.removal == 0 & sample.size == 0)
cort2.nuc2 <- subset(cort.m2, quality.removal == 0 & sample.size == 0)

ens_tubb <- subset(ens_tubb0,
                 quality.removal == 0 & sat.qual.test == 0 &
                 attempt.no != 4 & attempt.no != 14 & attempt.no != 6 & attempt.no != 19)

cort_tubb <- subset(cort_tubb0, qual.id == 0 & sat.qual == 0)
```


```{r ens area analysis,  warning = FALSE}
ens4$m.id <- paste(ens4$attempt,"-",ens4$well)
n.ens <- aggregate(data = ens4,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.ens <- aggregate(data = ens4,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.ens <- aggregate(data = ens4,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))

sum.avg.int.asyn.ens <- aggregate(data = ens4,
                          asyn_roi_tot_avg_int ~ Groups.f, FUN = mean)
sum.avg.int.pff.ens <- aggregate(data = ens4,
                          pff_roi_tot_avg_int ~ Groups.f, FUN = mean)

ens.replicates <- cbind(n.ens[,2],m.ens[,2],i.ens[,2],sum.avg.int.asyn.ens[,2],sum.avg.int.pff.ens[,2])
colnames(ens.replicates) <- c("N","m","no.images","asyn avg int","pff avg int")
rownames(ens.replicates) <-m.ens[,1]
ens.replicates

cort2$m.id <- paste(cort2$attempt,"-",cort2$well)
n.cort <- aggregate(data = cort2,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.cort <- aggregate(data = cort2,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.cort <- aggregate(data = cort2,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))

sum.avg.int.asyn.cort <- aggregate(data = cort2,
                          asyn_roi_tot_avg_int ~ Groups.f, FUN = mean)
sum.avg.int.pff.cort <- aggregate(data = cort2,
                          pff_roi_tot_avg_int ~ Groups.f, FUN = mean)

cort.replicates <- cbind(n.cort[,2],m.cort[,2],i.cort[,2],sum.avg.int.asyn.cort[,2],sum.avg.int.pff.cort[,2])
colnames(cort.replicates) <- c("N","m","no.images","asyn avg int","pff avg int")
rownames(cort.replicates) <-m.cort[,1]
cort.replicates

mod.ens.pff.p.area <- kruskal.test(pff_p_area_tot ~ Groups.f, data = ens4)
mod.ens.pff.p.area
ph.ens.pff.p.area <- pairwise.wilcox.test(ens4$pff_p_area_tot, ens4$Groups.f,p.adjust.method = "BH")
ph.ens.pff.p.area 

plot.ens.pff.area <- ggplot(ens4, aes(x=Groups.f, y=pff_p_area_tot, fill=Groups.f)) +
  ylab("Fraction of neural area with elevated PFF tag signal (ATTO-590)") + 
  xlab("Groups") + 
  labs(
    title = "Neuron PFF Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.5, tip_length = 0.05, vjust=0.4) +
    geom_signif(comparisons=list(c("2XPFF", "Control")), annotations="***",
              y_position = 0.65, tip_length = 0.05, vjust=0.4) +
    geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 0.75, tip_length = 0.05, vjust=0.4) +
    geom_signif(comparisons=list(c("6XPFF", "Control")), annotations="***",
              y_position = 0.95, tip_length = 0.05, vjust=0.4) +
    geom_signif(comparisons=list(c("PFF", "PFF+But")), annotations="0.076",
              y_position = 0.5, tip_length = 0.05, vjust=0) +
    geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="*",
              y_position = 0.58, tip_length = 0.05, vjust=0.4) +
    geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="**",
              y_position = 0.85, tip_length = 0.05, vjust=0.4) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_p_area_tot),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none")

plot.ens.pff.area

mod.ens.asyn.p.area <- kruskal.test(asyn_p_area_tot ~ Groups.f, data = ens4)
ph.ens.asyn.p.area <- pairwise.wilcox.test(ens4$asyn_p_area_tot, ens4$Groups.f,p.adjust.method = "BH")
ph.ens.asyn.p.area

plot.ens.asyn.area <- ggplot(ens4, aes(x=Groups.f, y=asyn_p_area_tot, fill=Groups.f)) +
  ylab("Fraction of neural area with elevated asyn") + 
  xlab("Groups") + 
  labs(
    title = "Neuron a-Syn Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 50, tip_length = 0.05, vjust=0.4) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_p_area_tot),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.001,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none")

plot.ens.asyn.area + coord_cartesian(ylim = c(0,0.01))
```


```{r ens area corrected for dapi,  warning = FALSE}
ens4.nuc$asyn_p_area_tot.corr <-ens4.nuc$asyn_p_area_tot/ens4.nuc$nuclei.no
ens4.nuc$pff_p_area_tot.corr <-ens4.nuc$pff_p_area_tot/ens4.nuc$nuclei.no

mod.ens.asyn.p.area.corr <- kruskal.test(asyn_p_area_tot.corr ~ Groups.f, data = ens4.nuc)
ph.ens.asyn.p.area.corr <- pairwise.wilcox.test(ens4.nuc$asyn_p_area_tot.corr, ens4.nuc$Groups.f,p.adjust.method = "BH")
ph.ens.asyn.p.area.corr 

mod.ens.pff.p.area.corr <- kruskal.test(pff_p_area_tot.corr ~ Groups.f, data = ens4.nuc)
ph.ens.pff.p.area.corr <- pairwise.wilcox.test(ens4.nuc$pff_p_area_tot.corr, ens4.nuc$Groups.f,p.adjust.method = "BH")
ph.ens.pff.p.area.corr 

plot.ens.asyn.area <- ggplot(ens4.nuc, aes(x=Groups.f, y=sqrt(asyn_p_area_tot.corr), fill=Groups.f)) +
  ylab("Sqrt Fraction Area per Nuclei (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "a-Syn Area per Nuclei (ENS)") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) +
  scale_fill_viridis_d()

plot.ens.asyn.area

plot.ens.pff.area <- ggplot(ens4.nuc, aes(x=Groups.f, y=sqrt(pff_p_area_tot.corr), fill=Groups.f)) +
  ylab("Sqrt Fraction Area per Nuclei (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "PFF Area per Nuclei (ENS)") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.05, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
    geom_signif(comparisons=list(c("2XPFF", "PFF")), annotations="***",
              y_position = 0.125, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
    geom_signif(comparisons=list(c("4XPFF", "2XPFF")), annotations="***",
              y_position = 0.20, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
    geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="**",
              y_position = 0.20, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) +
  coord_cartesian(ylim = c(0,0.28)) +
  scale_fill_viridis_d()

plot.ens.pff.area
    
mod.ens.asyn.nuc.int <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|attempt.no), 
                             data = ens4.nuc)

mod.ens.asyn.nuc.int2 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|well), 
                             data = ens4.nuc)

mod.ens.asyn.nuc.int3 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = ens4.nuc)

mod.ens.asyn.nuc.int4 <- lm(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no, 
                             data = ens4.nuc)

anova(mod.ens.asyn.nuc.int,mod.ens.asyn.nuc.int2,mod.ens.asyn.nuc.int3,mod.ens.asyn.nuc.int4)
anova(mod.ens.asyn.nuc.int,mod.ens.asyn.nuc.int3)
summary(mod.ens.asyn.nuc.int)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.asyn.avg.int.nuc <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + nuclei.no + 
        (1|attempt.no)
         , data = ens4.nuc))
summary(em.ens.asyn.avg.int.nuc <- emmeans(warp.asyn.avg.int.nuc, ~Groups.f, type = "response"))
pwpm(em.ens.asyn.avg.int.nuc)

em.ens.nuc.asyn.tot.int<- plot(em.ens.asyn.avg.int.nuc, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average Intensity adjusted by nuclei (a-Syn, ENS)",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.ens.nuc.asyn.tot.int

mod.ens.pff.nuc.int0 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no), 
                             data = ens4.nuc)

mod.ens.pff.nuc.int1 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|attempt.no), 
                             data = ens4.nuc)

mod.ens.pff.nuc.int2 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|well), 
                             data = ens4.nuc)

mod.ens.pff.nuc.int3 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = ens4.nuc)

mod.ens.pff.nuc.int4 <- lm(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + nuclei.no, 
                             data = ens4.nuc)

anova(mod.ens.pff.nuc.int0,mod.ens.pff.nuc.int1,mod.ens.pff.nuc.int2,mod.ens.pff.nuc.int3,mod.ens.pff.nuc.int4)
anova(mod.ens.pff.nuc.int1,mod.ens.pff.nuc.int3)
summary(mod.ens.pff.nuc.int1)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.pff.avg.int.nuc <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + nuclei.no + 
        (1|attempt.no)
         , data = ens4.nuc))
summary(em.ens.pff.avg.int.nuc <- emmeans(warp.pff.avg.int.nuc, ~Groups.f, type = "response"))
pwpm(em.ens.pff.avg.int.nuc)

em.ens.nuc.pff.tot.int<- plot(em.ens.pff.avg.int.nuc, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average Intensity adjusted by nuclei (PFF, ENS)",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.ens.nuc.pff.tot.int
```


```{r ens transformation and analysis avg int,  warning = FALSE}
#Identify whether transformation is needed
#skewness(ens4$asyn_roi_tot_avg_int, na.rm = TRUE)
hist(ens4$asyn_roi_tot_avg_int)
hist(ens4$t.asyn_roi_tot_avg_int <- 
       (ens4$asyn_roi_tot_avg_int)^(1/6))
skewness(ens4$t.asyn_roi_tot_avg_int, na.rm = TRUE)

#Analyze a-Syn
mod.ens.asyn.tot.int <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no), 
                             data = ens4)

mod.ens.asyn.tot.int2 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|well), 
                             data = ens4)

mod.ens.asyn.tot.int3 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = ens4)

mod.ens.asyn.tot.int4 <- lm(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f, 
                             data = ens4)

anova(mod.ens.asyn.tot.int,mod.ens.asyn.tot.int2,mod.ens.asyn.tot.int3,mod.ens.asyn.tot.int4)
anova(mod.ens.asyn.tot.int,mod.ens.asyn.tot.int3)
summary(mod.ens.asyn.tot.int3)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.asyn.avg.int.tot <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.ens.asyn.avg.int.tot <- emmeans(warp.asyn.avg.int.tot, ~Groups.f, type = "response"))
pwpm(em.ens.asyn.avg.int.tot)
plot(em.ens.asyn.avg.int.tot, comparisons=TRUE, type = "response")
#eff_size(em.ens.asyn.tot.int, sigma = sigma(mod.ens.asyn.tot.int), edf = 11)

#Analyze PFF
mod.ens.pff.tot.int <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no), 
                             data = ens4)

mod.ens.pff.tot.int2 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|well), 
                             data = ens4)

mod.ens.pff.tot.int3 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = ens4)

mod.ens.pff.tot.int4 <- lm(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f, 
                             data = ens4)

anova(mod.ens.pff.tot.int,mod.ens.pff.tot.int2,mod.ens.pff.tot.int3,mod.ens.pff.tot.int4)
anova(mod.ens.pff.tot.int,mod.ens.pff.tot.int3)
summary(mod.ens.pff.tot.int)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.pff.avg.int.tot <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens4))
summary(em.ens.pff.avg.int.tot <- emmeans(warp.pff.avg.int.tot, ~Groups.f, type = "response"))
pwpm(em.ens.pff.avg.int.tot)
plot(em.ens.pff.avg.int.tot, comparisons=TRUE, type = "response")
#eff_size(em.ens.pff.tot.int, sigma = sigma(mod.ens.pff.tot.int), edf = 11)

em.ens.asyn.tot.int<- plot(em.ens.asyn.avg.int.tot, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average a-Syn intensity within neuron (ENS)",
       x = "Estimated Marginal Mean (Average Intensity a-Syn)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.ens.asyn.tot.int

em.ens.pff.tot.int<- plot(em.ens.pff.avg.int.tot, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average PFF intensity within neuron (ENS)",
       x = "Estimated Marginal Mean (Average Intensity PFF)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.ens.pff.tot.int

# plot.ens.pff.int <- ggplot(ens4, aes(x=Groups.f, y=pff_roi_tot_avg_int, fill=Groups.f)) +
#   ylab("Average PFF Intensity within Neuron") + 
#   xlab("Groups") + 
#   labs(
#     title = "Neuron PFF Average Intensity") +
#   geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
#     stat_summary(fun.data = mean_se,  
#                  geom = "errorbar", position = 'dodge', colour = "black") + 
#   geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
#               y_position = 50, tip_length = 0.05, vjust=0.4) +
#   geom_signif(comparisons=list(c("2XPFF", "Control")), annotations="***",
#               y_position = 60, tip_length = 0.05, vjust=0.4) +
#   geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
#               y_position = 110, tip_length = 0.05, vjust=0) +
#   geom_signif(comparisons=list(c("6XPFF", "Control")), annotations="***",
#               y_position = 120, tip_length = 0.05, vjust=0) +
#   geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
#               y_position = 100, tip_length = 0.05, vjust=0) +
#   geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_roi_tot_avg_int),
#                alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
#                stackdir = "center", stackratio = 1) +
#   theme(legend.position="none")
# 
# plot.ens.pff.int
# 
# plot.ens.asyn.int <- ggplot(ens4, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=Groups.f)) +
#   ylab("Average a-Syn Intensity within Neuron") + 
#   xlab("Groups") + 
#   labs(
#     title = "Neuron a-Syn Average Intensity") +
#   geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
#     stat_summary(fun.data = mean_se,  
#                  geom = "errorbar", position = 'dodge', colour = "black") + 
#   geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_roi_tot_avg_int),
#                alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
#                stackdir = "center", stackratio = 1) +
#   geom_signif(comparisons=list(c("PFF", "Control")), annotations="**",
#               y_position = 40, tip_length = 0.05, vjust=0.4) +
#   geom_signif(comparisons=list(c("2XPFF", "Control")), annotations="**",
#               y_position = 45, tip_length = 0.05, vjust=0.4) +
#   geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
#               y_position = 50, tip_length = 0.05, vjust=0) +
#   geom_signif(comparisons=list(c("6XPFF", "Control")), annotations="***",
#               y_position = 55, tip_length = 0.05, vjust=0) +
#   geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="0.057",
#               y_position = 40, tip_length = 0.05, vjust=0) +
#   geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="**",
#               y_position = 60, tip_length = 0.05, vjust=0) +
#   theme(legend.position="none")
# 
# plot.ens.asyn.int

plot.ens.asyn.int.group <- ggplot(ens4, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=attempt.no)) +
  ylab("Average a-Syn Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Alpha synuclein average intensity ENS", 
    subtitle = "3-8 replicates, 8-21 wells, 40-90 images, Wilcox Non-parametric Test") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_roi_tot_avg_int,
                                groups = interaction(Groups.f, attempt.no)),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1)

plot.ens.asyn.int.group

ens3.simp <- subset(ens3, (Groups.f == "Control" | Groups.f == "PFF" | Groups.f == "But" | Groups.f == "PFF+But") )

plot.ens.pff.int.group <- ggplot(ens3.simp, aes(x=Groups.f, y=pff_roi_tot_avg_int, fill=attempt.no)) +
  ylab("Average PFF Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Variation exists within group and experimental replicate") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens3.simp, aes(x = Groups.f, y = pff_roi_tot_avg_int,
                                groups = interaction(Groups.f, attempt.no)),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1)

plot.ens.pff.int.group

plot.ens.asyn.int.group.simp <- ggplot(ens3, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=attempt.no)) +
  ylab("Average a-Syn Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Neuron a-Syn Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens3, aes(x = Groups.f, y = asyn_roi_tot_avg_int,
                                groups = interaction(Groups.f, attempt.no)),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1)

plot.ens.asyn.int.group.simp
```



```{r analysis by well data,  warning = FALSE}
#Setup
ens.w0 <- subset(ens0, 
                 quality.removal == 0 & sat.qual.test == 0 & pff.qual == 0 & but.qual == 0 &
                 attempt != 4 & attempt != 14 & attempt != 6)
                 
ens.w1 <- aggregate(cbind(
                    asyn_avg_int,
                    pff_avg_int,
                    asyn_roi_tot_avg_int,
                    pff_roi_tot_avg_int,
                    asyn_p_area_tot,
                    pff_p_area_tot) ~
                  group + attempt + well,
                  ens.w0, FUN = mean)

ens.w1$Groups.f <- factor(ens.w1$group, levels = c("Control","But", "LPS", "But+LPS",
                                             "PFF", "PFF+But", "PFF+LPS", "PFF+But+LPS",
                                             "2XPFF", "4XPFF", "4XPFF+But", "6XPFF"))
ens.w1$attempt.no <- factor(ens.w1$attempt)

#Avg Int
mod.ens.asyn.tot.int.w1 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f, data = ens.w1)
mod.ens.asyn.tot.int.w2 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + (1|attempt.no), data = ens.w1)
lrtest(mod.ens.asyn.tot.int.w1,mod.ens.asyn.tot.int.w2)
fixed.avg.int.tran <- make.tran("power", (1/6))
warp.asyn.avg.int.well <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens.w1))
summary(em.ens.asyn.avg.int.well <- emmeans(warp.asyn.avg.int.well, ~Groups.f, type = "response"))
pwpm(em.ens.asyn.avg.int.well)
plot(em.ens.asyn.avg.int.well, comparisons=TRUE, type = "response")

mod.ens.pff.tot.int.w1 <- lm(pff_roi_tot_avg_int^(1/6) ~ Groups.f, data = ens.w1)
mod.ens.pff.tot.int.w2 <- lmer(pff_roi_tot_avg_int^(1/6) ~ Groups.f + (1|attempt.no), data = ens.w1)
fixed.avg.int.tran <- make.tran("power", (1/6))
warp.pff.avg.int.well <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens.w1))
summary(em.ens.pff.avg.int.well <- emmeans(warp.pff.avg.int.well, ~Groups.f, type = "response"))
pwpm(em.ens.pff.avg.int.well)
plot(em.ens.pff.avg.int.well, comparisons=TRUE, type = "response")

#Non-parametric area
np.well.asyn.area <- pairwise.wilcox.test(ens.w1$asyn_p_area_tot, ens.w1$Groups.f,p.adjust.method = "BH")
np.well.asyn.area

np.well.pff.area <- pairwise.wilcox.test(ens.w1$pff_p_area_tot, ens.w1$Groups.f,p.adjust.method = "BH")
np.well.pff.area

plot.em.ens.well.asyn.int <- plot(em.ens.asyn.avg.int.well, 
                                 type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average a-Syn intensity averaged to well within neuron emmeans plot (ENS)",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups")
plot.em.ens.well.asyn.int

plot.em.ens.well.pff.int <- plot(em.ens.pff.avg.int.well, 
                                 type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average PFF tag intensity averaged to well within neuron emmeans plot (ENS)",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups")
plot.em.ens.well.pff.int

plot.ens.well.pff.int <- ggplot(ens.w1, aes(x=Groups.f, y=pff_roi_tot_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "PFF dosage results in an increase of a-Syn intensity", 
    subtitle = "3-8 replicates, 8-21 wells, 40-90 images, Wilcox Non-parametric Test") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens.w1, aes(x = Groups.f, y = pff_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none")

plot.ens.well.pff.int

plot.ens.well.asyn.int <- ggplot(ens.w1, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Alpha synuclein average intensity ENS", 
    subtitle = "3-8 replicates, 8-21 wells, 40-90 images, Wilcox Non-parametric Test") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens.w1, aes(x = Groups.f, y = asyn_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none")

plot.ens.well.asyn.int

```


```{r nuclei count by group ens,  warning = FALSE}
hist(ens4.nuc$nuclei.no)

mod.ens.nuclei.no <- kruskal.test(nuclei.no ~ Groups.f, data = ens4.nuc)
ph.ens.nuclei.no <- pairwise.wilcox.test(ens4.nuc$nuclei.no, ens4.nuc$Groups.f,p.adjust.method = "BH")
ph.ens.nuclei.no 

ens4.nuc$bin.nuclei <- cut(ens4.nuc$nuclei.no, breaks = c(0,5,15,30,60,120))
ens4.nuc$bin.nuclei <- factor(ens4.nuc$bin.nuclei)

mod.ens.nuclei.no2 <- lm(asyn_roi_tot_avg_int^(1/6) ~ bin.nuclei * Groups.f,data = ens4.nuc)
em.ens.nuclei.no2 <- emmeans(mod.ens.nuclei.no2, specs = "Groups.f", by="bin.nuclei")
pwpm(em.ens.nuclei.no2)

plot.ens.nuclei <- ggplot(ens4.nuc, aes(x=Groups.f, y=nuclei.no, fill=Groups.f)) +
  ylab("Number of Nuclei per image") + 
  xlab("Groups") + 
  labs(
    title = "Comparison of Nuclei Count per Groups")+
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens4.nuc, aes(x = Groups.f, y = nuclei.no),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none", text = element_text(size = 14))

plot.ens.nuclei


plot.ens.nuclei2 <- ggplot(ens4.nuc, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=bin.nuclei)) +
  ylab("Average Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Comparison of Nuclei Count per Groups")+
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens4.nuc, aes(x = Groups.f, y = asyn_roi_tot_avg_int,
                                groups = interaction(Groups.f, nuclei.no)),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(text = element_text(size = 14))

plot.ens.nuclei2

plot.ens.nuclei3 <- ggplot(ens4.nuc, aes(x=nuclei.no, y=asyn_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  xlab("Number of Nuclei in image") +
  ylab("Average a-Syn Intensity in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.nuclei3

ens4.nuc.2 <- subset(ens4.nuc, 
                 (Groups.f == "Control" | Groups.f == "But" | Groups.f == "PFF" | Groups.f == "PFF+But"))

ens4.nuc.3 <- subset(ens4.nuc, 
                 Groups.f == "Control")

plot.ens.nuclei3.2 <- ggplot(ens4.nuc.2, aes(x=nuclei.no, y=asyn_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. a-Syn intensity within Neuron (all groups)"
  ) +
  xlab("Number of Nuclei in image") +
  ylab("Average a-Syn Intensity in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.nuclei3.2

plot.ens.nuclei3.2p <- ggplot(ens4.nuc.2, aes(x=nuclei.no, y=pff_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. PFF intensity within Neuron (all groups)",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Average PFF Intensity in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.nuclei3.2p

plot.ens.comp3.3 <- ggplot(ens4.nuc.2, aes(x=pff_roi_tot_avg_int, y=asyn_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Comparison of a-Syn intensity with PFF intensity",
    col = "Groups") +
  xlab("Average PFF Intensity in Neuron") +
  ylab("Average a-Syn Intensity in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.comp3.3

#Scatterplot Sample Size
ens4.nuc.2$m.id <- paste(ens4.nuc.2$attempt,"-",ens4.nuc.2$well)
n.ens.nuc <- aggregate(data = ens4.nuc.2,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.ens.nuc <- aggregate(data = ens4.nuc.2,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.ens.nuc <- aggregate(data = ens4.nuc.2,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))
ens.nuc.replicates <- cbind(n.ens.nuc[,2],m.ens.nuc[,2],i.ens.nuc[,2])
colnames(ens.nuc.replicates) <- c("N","m","no.images")
rownames(ens.nuc.replicates) <- c("Control","But","PFF","PFF+But")
ens.nuc.replicates

#Corrected Plots Sample Size
ens4.nuc$m.id <- paste(ens4.nuc$attempt,"-",ens4.nuc$well)
n.ens.nuc <- aggregate(data = ens4.nuc,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.ens.nuc <- aggregate(data = ens4.nuc,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.ens.nuc <- aggregate(data = ens4.nuc,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))
ens.nuc.replicates.corr <- cbind(n.ens.nuc[,2],m.ens.nuc[,2],i.ens.nuc[,2])
colnames(ens.nuc.replicates.corr) <- c("N","m","no.images")
rownames(ens.nuc.replicates.corr) <- n.ens.nuc[,1]
ens.nuc.replicates.corr

plot.ens.nuclei3.3p <- ggplot(ens4.nuc.3, aes(x=nuclei.no, y=pff_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. PFF intensity within Neuron",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Average PFF Intensity in Neuron") +
  theme(legend.position="none", text = element_text(size = 14))

plot.ens.nuclei3.3p

plot.ens.nuclei3.3a <- ggplot(ens4.nuc.3, aes(x=nuclei.no, y=asyn_roi_tot_avg_int, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. a-Syn intensity within Neuron",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Average a-Syn Intensity within Neuron") +
  theme(legend.position="none", text = element_text(size = 14))

plot.ens.nuclei3.3a

plot.ens.nuclei3.4p <- ggplot(ens4.nuc.3, aes(x=nuclei.no, y=pff_p_area_tot, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. PFF elevated area",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Fraction of neural area with elevated PFF (ATTO-590)") +
  theme(legend.position="none",text = element_text(size = 14))

plot.ens.nuclei3.4p

plot.ens.nuclei3.4a <- ggplot(ens4.nuc.3, aes(x=nuclei.no, y=asyn_p_area_tot, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. a-Syn elevated area",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Fraction of neural area with elevated a-Syn") +
  theme(legend.position="none",text = element_text(size = 14))

plot.ens.nuclei3.4a

plot.ens.nuclei3.5p <- ggplot(ens4.nuc.2, aes(x=nuclei.no, y=pff_p_area_tot, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. fraction elevated PFF area with Neuron (all groups)",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Average PFF elevated area in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.nuclei3.5p

plot.ens.nuclei3.5a <- ggplot(ens4.nuc.2, aes(x=nuclei.no, y=asyn_p_area_tot, colour=Groups.f)) +
  geom_point() +
  labs(
    title = "Nuclei count vs. fraction elevated a-Syn area with Neuron (all groups)",
    col = "Groups") +
  xlab("Number of Nuclei in image") +
  ylab("Average a-Syn elevated area in Neuron") +
  theme(text = element_text(size = 14))

plot.ens.nuclei3.5a


mod.ens.nuclei.no3a <- lm(asyn_roi_tot_avg_int^(1/6) ~ nuclei.no,data = ens4.nuc.3)
summary(mod.ens.nuclei.no3a)


mod.ens.nuclei.no3p <- lm(pff_roi_tot_avg_int^(1/6) ~ nuclei.no,data = ens4.nuc.3)
summary(mod.ens.nuclei.no3p)

mod.ens.nuclei.no2b <- lm(asyn_roi_tot_avg_int^(1/6) ~ bin.nuclei * Groups.f,data = ens4.nuc.2)
em.ens.nuclei.no2b <- emmeans(mod.ens.nuclei.no2b, specs = "Groups.f", by="bin.nuclei")
pwpm(em.ens.nuclei.no2b)


mod.ens.nuc.asyn.tot.int <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no,data = ens4.nuc)
summary(mod.ens.nuclei.no3a)


mod.ens.nuclei.no3p <- lm(pff_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no,data = ens4.nuc)
summary(mod.ens.nuclei.no3p)
```

```{r nuclei adjusted,  warning = FALSE}
#Subset
ens5.nuc <- subset(ens4.nuc, nuclei.no < 20)

#Avg Int
mod.asyn.ens.nuc.int <- lmer(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens5.nuc)
warp.ens.asyn.avg.int.raw <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_line_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens5.nuc))
summary(em.asyn.ens.nuc.int <- emmeans(warp.ens.asyn.avg.int.raw, ~Groups.f, type = "response"))
pwpm(em.asyn.ens.nuc.int)

mod.pff.ens.nuc.int <- lmer(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens5.nuc)
warp.ens.pff.avg.int.raw <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_line_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens5.nuc))
summary(em.pff.ens.nuc.int <- emmeans(warp.ens.pff.avg.int.raw, ~Groups.f, type = "response"))
pwpm(em.pff.ens.nuc.int)

#Area
mod.ens5.nuc.asyn.p.area.corr <- kruskal.test(asyn_p_area_tot.corr ~ Groups.f, data = ens5.nuc)
ph.ens5.nuc.asyn.p.area.corr <- pairwise.wilcox.test(ens5.nuc$asyn_p_area_tot.corr, ens5.nuc$Groups.f,p.adjust.method = "BH")
ph.ens5.nuc.asyn.p.area.corr 

mod.ens5.nuc.pff.p.area <- kruskal.test(pff_p_area_tot ~ Groups.f, data = ens5.nuc)
ph.ens5.nuc.pff.p.area <- pairwise.wilcox.test(ens5.nuc$pff_p_area_tot, ens5.nuc$Groups.f,p.adjust.method = "BH")
ph.ens5.nuc.pff.p.area 

```


```{r cortical,  warning = FALSE}
#Sample Size
cort2$m.id <- paste(cort2$attempt,"-",cort2$well)
n.cort <- aggregate(data = cort2,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.cort <- aggregate(data = cort2,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.cort <- aggregate(data = cort2,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))
cort.replicates <- cbind(n.cort[,2],m.cort[,2],i.cort[,2])
colnames(cort.replicates) <- c("N","m","no.images")
rownames(cort.replicates) <- c("Control","PFF","4XPFF")
cort.replicates

#Transformation
hist(cort2$pff_roi_tot_avg_int)
hist(cort2$t.pff_roi_tot_avg_int <-
       (cort2$pff_roi_tot_avg_int)^(1/6))

#Analyze a-Syn
mod.cor.asyn.tot.int <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no), 
                             data = cort2)

mod.cor.asyn.tot.int2 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|well), 
                             data = cort2)

mod.cor.asyn.tot.int3 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = cort2)

mod.cor.asyn.tot.int4 <- lm(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f, 
                             data = cort2)

anova(mod.cor.asyn.tot.int,mod.cor.asyn.tot.int2,mod.cor.asyn.tot.int3,mod.cor.asyn.tot.int4)
anova(mod.cor.asyn.tot.int,mod.cor.asyn.tot.int3)
summary(mod.cor.asyn.tot.int)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.cort.asyn.avg.int.tot <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = cort2))
summary(em.cor.asyn.avg.int.tot <- emmeans(warp.cort.asyn.avg.int.tot, ~Groups.f, type = "response"))
pwpm(em.cor.asyn.avg.int.tot)
plot(em.cor.asyn.avg.int.tot, comparisons=TRUE, type = "response")
#eff_size(em.cor.asyn.avg.int.tot, sigma = sigma(mod.cor.asyn.tot.int), edf = 2)

#Analyze PFF
mod.cor.pff.tot.int <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no), 
                             data = cort2)

mod.cor.pff.tot.int2 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|well), 
                             data = cort2)

mod.cor.pff.tot.int3 <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                               (1|attempt.no) + (1|attempt.no:well), 
                             data = cort2)

mod.cor.pff.tot.int4 <- lm(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f, 
                             data = cort2)

anova(mod.cor.pff.tot.int,mod.cor.pff.tot.int2,mod.cor.pff.tot.int3,mod.cor.pff.tot.int4)
anova(mod.cor.pff.tot.int,mod.cor.pff.tot.int3)
summary(mod.cor.pff.tot.int)

fixed.avg.int.tran <- make.tran("power", (1/6))
warp.pff.avg.int.tot <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = cort2))
summary(em.cor.pff.avg.int.tot <- emmeans(warp.pff.avg.int.tot, ~Groups.f, type = "response"))
pwpm(em.cor.pff.avg.int.tot)
plot(em.cor.pff.avg.int.tot, comparisons=TRUE, type = "response")
#eff_size(em.cor.pff.tot.int, sigma = sigma(mod.cor.pff.tot.int), edf = 2)

em.cor.asyn.tot.int<- plot(em.cor.asyn.avg.int.tot, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average a-Syn intensity within neuron (Cortex)",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.cor.asyn.tot.int

em.cor.pff.tot.int<- plot(em.cor.pff.avg.int.tot, type = "response", comparisons = TRUE, plotIT=FALSE) +
  labs(title = "Average PFF tag intensity within neuron (Cortex)",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.cor.pff.tot.int

plot.cort.asyn.int <- ggplot(cort2, aes(x=Groups.f, y=asyn_roi_tot_avg_int)) +
  ylab("Average a-Syn Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Average a-Syn intensity within neuron (Cortex)") +
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean',fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 75, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 82, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="**",
              y_position = 75, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = cort2, aes(x = Groups.f, y = asyn_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.6,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank())

plot.cort.asyn.int

plot.cort.pff.int <- ggplot(cort2, aes(x=Groups.f, y=pff_roi_tot_avg_int)) +
  ylab("Average PFF (ATTO-590) Intensity within Neuron") + 
  xlab("Groups") + 
  labs(
    title = "Cortical Average Intensity with Neuron (PFF)") +
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean',fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 80, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 87, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
              y_position = 80, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = cort2, aes(x = Groups.f, y = pff_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.35,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank())

plot.cort.pff.int
```

```{r cort area,  warning = FALSE}
hist(cort2$asyn_p_area_tot)
hist(cort2$t.asyn_p_area_tot <-
       1/(cort2$asyn_p_area_tot+0.00001))
skewness(cort2$t.asyn_p_area_tot, na.rm = FALSE)

mod.cor.asyn.p.area <- kruskal.test(asyn_p_area_tot ~ Groups.f, data = cort2)
ph.cor.asyn.p.area <- pairwise.wilcox.test(cort2$asyn_p_area_tot, cort2$Groups.f,p.adjust.method = "BH")
ph.cor.asyn.p.area 

mod.cor.pff.p.area <- kruskal.test(pff_p_area_tot ~ Groups.f, data = cort2)
ph.cor.pff.p.area <- pairwise.wilcox.test(cort2$pff_p_area_tot, cort2$Groups.f,p.adjust.method = "BH")
ph.cor.pff.p.area 

plot.cort.asyn.area <- ggplot(cort2, aes(x=Groups.f, y=sqrt(asyn_p_area_tot))) +
  ylab("Sqrt Fraction Area within Neuron (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "Cortex a-Syn Area")+
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean', fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black")  +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 0.2, tip_length = 0.01, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
              y_position = 0.23, tip_length = 0.01, size = 0.8, textsize = 8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank())

plot.cort.asyn.area 

# cort3 <- subset(cort2, attempt.no != 17)
# plot.cort.asyn.area <- ggplot(cort3, aes(x=Groups.f, y=asyn_p_area_tot, fill=Groups.f)) +
#   scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
#   geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
#     stat_summary(fun.data = mean_se, geom = "errorbar", position = 'dodge', colour = "black")
# 
# plot.cort.asyn.area

plot.cort.pff.area <- ggplot(cort2, aes(x=Groups.f, y=sqrt(pff_p_area_tot))) +
  ylab("Sqrt Fraction Area within Neuron (PFF)") + 
  xlab("Groups") + 
  labs(
    title = "Cortex PFF Area")+
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean', fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.625, tip_length = 0.05, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 0.725, tip_length = 0.05, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
              y_position = 0.625, tip_length = 0.05, size = 0.8, textsize = 8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank()) 

plot.cort.pff.area

```

```{r cort avg int corrected for nuclei}
#a-Syn
mod.cor.asyn.tot.int.corr <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no, data = cort2.nuc)
mod.cor.asyn.tot.int.corr2 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|attempt.no),
                                   data = cort2.nuc)
mod.cor.asyn.tot.int.corr3 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|well),
                                   data = cort2.nuc)
mod.cor.asyn.tot.int.corr4 <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|attempt.no) + (1|attempt:well),
                                   data = cort2.nuc)
lrtest(mod.cor.asyn.tot.int.corr,mod.cor.asyn.tot.int.corr2,
      mod.cor.asyn.tot.int.corr3,mod.cor.asyn.tot.int.corr4)
anova(mod.cor.asyn.tot.int.corr2,mod.cor.asyn.tot.int.corr3,mod.cor.asyn.tot.int.corr4)

warp.cor.asyn.tot.int.corr <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + nuclei.no +
        (1|attempt.no) + (1|attempt:well)
         , data = cort2.nuc))
summary(cor.asyn.tot.int.corr <- emmeans(warp.cor.asyn.tot.int.corr, ~Groups.f, type = "response"))
pwpm(cor.asyn.tot.int.corr)

em.cor.asyn.tot.int.corr <- emmeans(mod.cor.asyn.tot.int.corr, ~"Groups.f", )
pwpm(em.cor.asyn.tot.int.corr)

em.cor.asyn.tot.int.corr<- plot(em.cor.asyn.tot.int.corr, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Average Intensity adjusted by nuclei (a-Syn, Cortex)",
       x = "Estimated Marginal Mean (Average Intensity a-Syn)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.cor.asyn.tot.int.corr

#PFF
mod.cor.pff.tot.int.corr <- lm(pff_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no, data = cort2.nuc)
mod.cor.pff.tot.int.corr2 <- lmer(pff_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|attempt.no),
                                   data = cort2.nuc)
mod.cor.pff.tot.int.corr3 <- lmer(pff_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|well),
                                   data = cort2.nuc)
mod.cor.pff.tot.int.corr4 <- lmer(pff_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no +
                                     (1|attempt.no) + (1|attempt:well),
                                   data = cort2.nuc)
lrtest(mod.cor.pff.tot.int.corr,mod.cor.pff.tot.int.corr2,
      mod.cor.pff.tot.int.corr3,mod.cor.pff.tot.int.corr4)
anova(mod.cor.pff.tot.int.corr2,mod.cor.pff.tot.int.corr3,mod.cor.pff.tot.int.corr4)

warp.cor.pff.tot.int.corr <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + nuclei.no +
        (1|attempt.no) + (1|attempt:well)
         , data = cort2.nuc))
summary(cor.pff.tot.int.corr <- emmeans(warp.cor.pff.tot.int.corr, ~Groups.f, type = "response"))
pwpm(cor.pff.tot.int.corr)

em.cor.pff.tot.int.corr <- emmeans(mod.cor.pff.tot.int.corr, ~"Groups.f", )
pwpm(em.cor.pff.tot.int.corr)

em.cor.pff.tot.int.corr<- plot(em.cor.pff.tot.int.corr, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Average Intensity adjusted by nuclei (PFF, Cortex)",
       x = "Estimated Marginal Mean (Average Intensity PFF)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())
em.cor.pff.tot.int.corr

```

```{r cort area corrected for nuclei,  warning = FALSE}
cort2.nuc$asyn_p_area_tot.corr <-cort2.nuc$asyn_p_area_tot/cort2.nuc$nuclei.no
cort2.nuc$pff_p_area_tot.corr <-cort2.nuc$pff_p_area_tot/cort2.nuc$nuclei.no

mod.cor.asyn.p.area.corr <- kruskal.test(asyn_p_area_tot.corr ~ Groups.f, data = cort2.nuc)
ph.cor.asyn.p.area.corr <- pairwise.wilcox.test(
    cort2.nuc$asyn_p_area_tot.corr, 
    cort2.nuc$Groups.f,
    p.adjust.method = "BH")
ph.cor.asyn.p.area.corr 

mod.cor.pff.p.area.corr <- kruskal.test(pff_p_area_tot.corr ~ Groups.f, data = cort2.nuc)
ph.cor.pff.p.area.corr <- pairwise.wilcox.test(
    cort2.nuc$pff_p_area_tot.corr, 
    cort2.nuc$Groups.f,
    p.adjust.method = "BH")
ph.cor.pff.p.area.corr 

plot.cort.asyn.area <- ggplot(cort2.nuc, aes(x=Groups.f, y=sqrt(asyn_p_area_tot.corr))) +
  ylab("Sqrt Fraction Area with Neuron per nuclei (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "a-Syn Area per Nuclei (Cortex)")+
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean',fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="0.06",
              y_position = 0.04, tip_length = 0.005, size = 0.8, textsize = 4, vjust=0) +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 0.05, tip_length = 0.005, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
              y_position = 0.04, tip_length = 0.005, size = 0.8, textsize = 8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.09))

plot.cort.asyn.area

# cort3 <- subset(cort2.nuc, attempt.no != 17)
# plot.cort.asyn.area <- ggplot(cort3, aes(x=Groups.f, y=asyn_p_area_tot.corr, fill=Groups.f)) +
#   scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
#   geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
#     stat_summary(fun.data = mean_se, geom = "errorbar", position = 'dodge', colour = "black")
# 
# plot.cort.asyn.area

plot.cort.pff.area <- ggplot(cort2.nuc, aes(x=Groups.f, y=sqrt(pff_p_area_tot.corr))) +
  ylab("Sqrt Fraction Area with Neuron per nuclei (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "PFF Area per Nuclei (Cortex)")+
  scale_x_discrete(labels=c("Control", "1 µg PFF", "4µg PFF")) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean',fill = "aquamarine4") +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.30, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
              y_position = 0.35, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
              y_position = 0.30, tip_length = 0.02, size = 0.8, textsize = 8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.title.x = element_blank()) 
  
plot.cort.pff.area

cort2.nuc$m.id <- paste(cort2.nuc$attempt,"-",cort2.nuc$well)
n.cort.nuc <- aggregate(data = cort2.nuc,
                          attempt ~ Groups.f,
                          function(attempt) length(unique(attempt)))
m.cort.nuc <- aggregate(data = cort2.nuc,
                          m.id ~ Groups.f,
                          function(m.id) length(unique(m.id)))
i.cort.nuc <- aggregate(data = cort2.nuc,
                          number.analyzed ~ Groups.f,
                          function(number.analyzed) length(unique(number.analyzed)))
cort.nuc.replicates.corr <- cbind(n.cort.nuc[,2],m.cort.nuc[,2],i.cort.nuc[,2])
colnames(cort.nuc.replicates.corr) <- c("N","m","no.images")
rownames(cort.nuc.replicates.corr) <- n.cort.nuc[,1]
cort.nuc.replicates.corr


```


```{r ens optimization raw}
#Area
mod.asyn.p.raw.area <- kruskal.test(asyn_p_area ~ Groups.f, data = ens4)
ph.asyn.p.raw.area <- pairwise.wilcox.test(ens4$asyn_p_area, ens4$Groups.f,p.adjust.method = "BH")
ph.asyn.p.raw.area

mod.pff.p.raw.area  <- kruskal.test(pff_p_area ~ Groups.f, data = ens4)
ph.pff.p.raw.area <- pairwise.wilcox.test(ens4$pff_p_area, ens4$Groups.f,p.adjust.method = "BH")
ph.pff.p.raw.area

plot.ens.asyn.raw.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(asyn_p_area), fill=Groups.f)) +
  ylab("Sqrt fraction area within whole image (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "Whole Image a-Syn Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("Control", "4XPFF")), annotations="*",
              y_position = 0.0035, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="*",
              y_position = 0.006, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="*",
              y_position = 0.001, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.03)) +
  scale_fill_viridis_d()

plot.ens.asyn.raw.area

plot.ens.pff.raw.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(pff_p_area), fill=Groups.f)) +
  ylab("Sqrt fraction area within whole image (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "Whole Image PFF Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.1, tip_length = 0.01, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "PFF")), annotations="***",
              y_position = 0.15, tip_length = 0.01, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 0.24, tip_length = 0.01, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 0.1, tip_length = 0.01, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 0.24, tip_length = 0.01, size = 0.8, textsize=8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  scale_fill_viridis_d()

plot.ens.pff.raw.area

#Avg Int Raw
mod.asyn.raw.int0 <- lm(asyn_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.asyn.raw.int1 <- lmer(asyn_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.asyn.raw.int2 <- lmer(asyn_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.asyn.raw.int3 <- lmer(asyn_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.asyn.raw.int0,mod.asyn.raw.int1,mod.asyn.raw.int2,mod.asyn.raw.int3)
anova(mod.asyn.raw.int1,mod.asyn.raw.int2,mod.asyn.raw.int3)

mod.asyn.raw.int <- lmer(asyn_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens4)
warp.ens.asyn.avg.int.raw <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.asyn.raw.int <- emmeans(warp.ens.asyn.avg.int.raw, ~Groups.f, type = "response"))
pwpm(em.asyn.raw.int)

mod.pff.raw.int0 <- lm(pff_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.pff.raw.int1 <- lmer(pff_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.pff.raw.int2 <- lmer(pff_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.pff.raw.int3 <- lmer(pff_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.pff.raw.int0,mod.pff.raw.int1,mod.pff.raw.int2,mod.pff.raw.int3)
anova(mod.pff.raw.int1,mod.pff.raw.int2,mod.pff.raw.int3)

mod.pff.raw.int <- lmer(pff_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens4)
warp.ens.pff.avg.int.raw <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.pff.raw.int <- emmeans(warp.ens.pff.avg.int.raw, ~Groups.f, type = "response"))
pwpm(em.pff.raw.int)

plot.ens.asyn.raw.int <- ggplot(ens4, aes(x=Groups.f, y=asyn_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Whole Image a-Syn Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("LPS", "Control")), annotations="***",
              y_position = 28, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("LPS", "But+LPS")), annotations="**",
              y_position = 28, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "PFF")), annotations="***",
              y_position = 32, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "2XPFF")), annotations="***",
              y_position = 36, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="***",
              y_position = 40, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  coord_cartesian(ylim = c(0,43)) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  scale_fill_viridis_d()

plot.ens.asyn.raw.int

plot.ens.pff.raw.int <- ggplot(ens4, aes(x=Groups.f, y=pff_avg_int, fill=Groups.f)) +
  ylab("Average PFF (ATTO-590) Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Whole Image PFF Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 20, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 20, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 35, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "Control")), annotations="***",
              y_position = 28, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 35, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "6XPFF")), annotations="**",
              y_position = 39, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,42)) +
  scale_fill_viridis_d()

plot.ens.pff.raw.int

em.ens.asyn.raw.int<- plot(em.asyn.raw.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Whole Image a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups")
em.ens.asyn.raw.int

em.ens.pff.raw.int<- plot(em.pff.raw.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Whole Image PFF Average Intensity",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups")
em.ens.pff.raw.int
```

```{r ens optimization neurite}
#Area
mod.asyn.p.line.area <- kruskal.test(asyn_p_area_line ~ Groups.f, data = ens4)
ph.asyn.p.line.area <- pairwise.wilcox.test(ens4$asyn_p_area_line, ens4$Groups.f,p.adjust.method = "BH")
ph.asyn.p.line.area

mod.pff.p.line.area  <- kruskal.test(pff_p_area_line ~ Groups.f, data = ens4)
ph.pff.p.line.area <- pairwise.wilcox.test(ens4$pff_p_area_line, ens4$Groups.f,p.adjust.method = "BH")
ph.pff.p.line.area

plot.ens.asyn.line.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(asyn_p_area_line), fill=Groups.f)) +
  ylab("Sqrt fraction area with neurite (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "Neurite a-Syn Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
  geom_signif(comparisons=list(c("Control", "But")), annotations="*",
              y_position = 0.001, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "4XPFF")), annotations="*",
              y_position = 0.01, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="*",
              y_position = 0.013, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.038)) +
  scale_fill_viridis_d()

plot.ens.asyn.line.area

plot.ens.pff.line.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(pff_p_area_line), fill=Groups.f)) +
  ylab("Sqrt fraction area with neurite (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "Neurite PFF Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.2, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "PFF")), annotations="*",
              y_position = 0.3, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 0.48, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="**",
              y_position = 0.48, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But")), annotations="0.0546",
              y_position = 0.15, tip_length = 0.015, size = 0.8, textsize=4, vjust=0) +
  geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="*",
              y_position = 0.2, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="**",
              y_position = 0.25, tip_length = 0.015, size = 0.8, textsize=8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.55)) +
  scale_fill_viridis_d()
  
plot.ens.pff.line.area

#Avg Int Neurite
mod.asyn.line.int0 <- lm(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.asyn.line.int1 <- lmer(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.asyn.line.int2 <- lmer(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.asyn.line.int3 <- lmer(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.asyn.line.int0,mod.asyn.line.int1,mod.asyn.line.int2,mod.asyn.line.int3)
anova(mod.asyn.line.int1,mod.asyn.line.int2,mod.asyn.line.int3)

mod.asyn.line.int <- lmer(asyn_roi_line_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens4)
warp.ens.asyn.avg.int.line <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_line_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.asyn.line.int <- emmeans(warp.ens.asyn.avg.int.line, ~Groups.f, type = "response"))
pwpm(em.asyn.line.int)

mod.pff.line.int0 <- lm(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.pff.line.int1 <- lmer(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.pff.line.int2 <- lmer(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.pff.line.int3 <- lmer(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.pff.line.int0,mod.pff.line.int1,mod.pff.line.int2,mod.pff.line.int3)
anova(mod.pff.line.int1,mod.pff.line.int2,mod.pff.line.int3)

mod.pff.line.int <- lmer(pff_roi_line_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens4)
warp.ens.pff.avg.int.line <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_line_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.pff.line.int <- emmeans(warp.ens.pff.avg.int.line, ~Groups.f, type = "response"))
pwpm(em.pff.line.int)

plot.ens.asyn.neurite.int <- ggplot(ens4, aes(x=Groups.f, y=asyn_roi_line_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Neurite a-Syn Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("LPS", "Control")), annotations="**",
              y_position = 50, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("LPS", "But+LPS")), annotations="**",
              y_position = 50, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="*",
              y_position = 60, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFFl", "PFF+But")), annotations="*",
              y_position = 50, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_roi_line_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  scale_fill_viridis_d()

plot.ens.asyn.neurite.int

plot.ens.pff.neurite.int <- ggplot(ens4, aes(x=Groups.f, y=pff_roi_line_avg_int, fill=Groups.f)) +
  ylab("Average PFF (ATTO-590) Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Neurite PFF Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 60, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "2XPFF")), annotations="***",
              y_position = 75, tip_length = 0.02, size = 0.8, textsize=8, vjust=0) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 120, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 120, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "6XPFF")), annotations="***",
              y_position = 135, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_roi_line_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,145)) +
  scale_fill_viridis_d()

plot.ens.pff.neurite.int

em.ens.asyn.line.int<- plot(em.asyn.line.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Neurite a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups")
em.ens.asyn.line.int

em.ens.pff.line.int<- plot(em.pff.line.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Neurite a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups")
em.ens.pff.line.int
```

```{r ens optimization soma}
#Area
mod.asyn.p.soma.area <- kruskal.test(asyn_p_area_soma ~ Groups.f, data = ens4)
ph.asyn.p.soma.area <- pairwise.wilcox.test(ens4$asyn_p_area_soma, ens4$Groups.f,p.adjust.method = "BH")
ph.asyn.p.soma.area

mod.pff.p.soma.area  <- kruskal.test(pff_p_area_soma ~ Groups.f, data = ens4)
ph.pff.p.soma.area <- pairwise.wilcox.test(ens4$pff_p_area_soma, ens4$Groups.f,p.adjust.method = "BH")
ph.pff.p.soma.area

plot.ens.asyn.soma.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(asyn_p_area_soma), fill=Groups.f)) +
  ylab("Sqrt fraction area within soma (a-Syn)") + 
  xlab("Groups") + 
  labs(
    title = "Soma a-Syn Area") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
  geom_signif(comparisons=list(c("Control", "4XPFF")), annotations="**",
              y_position = 0.025, tip_length = 0.0025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="*",
              y_position = 0.015, tip_length = 0.0025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="***",
              y_position = 0.035, tip_length = 0.0025, size = 0.8, textsize=8, vjust=0.4) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.08)) +
  scale_fill_viridis_d()

plot.ens.asyn.soma.area

plot.ens.pff.soma.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(pff_p_area_soma), fill=Groups.f)) +
  ylab("Sqrt fraction area within soma (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "Soma PFF Area") +  
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.27, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "PFF")), annotations="***",
              y_position = 0.43, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 0.65, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="***",
              y_position = 0.27, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 0.35, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 0.65, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("6XPFF", "4XPFF")), annotations="***",
              y_position = 0.75, tip_length = 0.025, size = 0.8, textsize=8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.85)) +
  scale_fill_viridis_d()

plot.ens.pff.soma.area

#Avg Int Soma
mod.asyn.soma.int0 <- lm(asyn_roi_soma_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.asyn.soma.int1 <- lmer(asyn_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.asyn.soma.int2 <- lmer(asyn_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.asyn.soma.int3 <- lmer(asyn_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.asyn.soma.int0,mod.asyn.soma.int1,mod.asyn.soma.int2,mod.asyn.soma.int3)
anova(mod.asyn.soma.int1,mod.asyn.soma.int2,mod.asyn.soma.int3)

warp.ens.asyn.avg.int.soma <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_soma_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.asyn.soma.int <- emmeans(warp.ens.asyn.avg.int.soma, ~Groups.f, type = "response"))
pwpm(em.asyn.soma.int)

mod.pff.soma.int0 <- lm(pff_roi_soma_avg_int^(1/6) ~ 1 + Groups.f, data = ens4)
mod.pff.soma.int1 <- lmer(pff_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no), data = ens4)
mod.pff.soma.int2 <- lmer(pff_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|well), data = ens4)
mod.pff.soma.int3 <- lmer(pff_roi_soma_avg_int^(1/6) ~ 1 + Groups.f + 
                            (1|attempt.no) + (1|attempt.no:well), data = ens4)

lrtest(mod.pff.soma.int0,mod.pff.soma.int1,mod.pff.soma.int2,mod.pff.soma.int3)
anova(mod.pff.soma.int1,mod.pff.soma.int2,mod.pff.soma.int3)

warp.ens.pff.avg.int.soma <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_soma_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.pff.soma.int <- emmeans(warp.ens.pff.avg.int.soma, ~Groups.f, type = "response"))
pwpm(em.pff.soma.int)

plot.ens.asyn.soma.int <- ggplot(ens4, aes(x=Groups.f, y=asyn_roi_soma_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Soma a-Syn Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 100, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_roi_soma_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  scale_fill_viridis_d()

plot.ens.asyn.soma.int

plot.ens.pff.soma.int <- ggplot(ens4, aes(x=Groups.f, y=pff_roi_soma_avg_int, fill=Groups.f)) +
  ylab("Average PFF (ATTO-590) Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Soma PFF Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 100, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="*",
              y_position = 100, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "2XPFF")), annotations="***",
              y_position = 115, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 135, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 135, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("6XPFF", "4XPFF")), annotations="***",
              y_position = 150, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 135, tip_length = 0.03, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_roi_soma_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,160)) +
  scale_fill_viridis_d()

plot.ens.pff.soma.int

em.ens.asyn.soma.int<- plot(em.asyn.soma.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Soma a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups")
em.ens.asyn.soma.int

em.ens.pff.soma.int<- plot(em.pff.soma.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Soma a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity PFF)",
       y = "Groups")
em.ens.pff.soma.int
```


```{r ens optimization neuron}
#Area
mod.asyn.p.neuron.area <- kruskal.test(asyn_p_area_tot ~ Groups.f, data = ens4)
ph.asyn.p.neuron.area <- pairwise.wilcox.test(ens4$asyn_p_area_tot, ens4$Groups.f,p.adjust.method = "BH")
ph.asyn.p.neuron.area

mod.pff.p.neuron.area  <- kruskal.test(pff_p_area_tot ~ Groups.f, data = ens4)
ph.pff.p.neuron.area <- pairwise.wilcox.test(ens4$pff_p_area_tot, ens4$Groups.f,p.adjust.method = "BH")
ph.pff.p.neuron.area

plot.ens.asyn.neuron.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(asyn_p_area_tot), fill=Groups.f)) +
  ylab("Sqrt fraction area within neuron (a-Syn") + 
  xlab("Groups") + 
  labs(
    title = "Neuron a-Syn Area") +
  coord_cartesian(ylim = c(0,0.05)) +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
  geom_signif(comparisons=list(c("Control", "But")), annotations="*",
              y_position = 0.001, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "4XPFF")), annotations="*",
              y_position = 0.008, tip_length = 0.002, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "6XPFF")), annotations="0.0522",
              y_position = 0.015, tip_length = 0.002, size = 0.8, textsize=4, vjust=0) +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) +
  scale_fill_viridis_d()

plot.ens.asyn.neuron.area

plot.ens.pff.neuron.area <- ggplot(ens4, aes(x=Groups.f, y=sqrt(pff_p_area_tot), fill=Groups.f)) +
  ylab("Sqrt fraction area within neuron (PFF tag)") + 
  xlab("Groups") + 
  labs(
    title = "Neuron PFF Area") +  
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "BLACK") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 0.2, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But")), annotations="0.0581",
              y_position = 0.2, tip_length = 0.02, size = 0.8, textsize=4, vjust=0) +
  geom_signif(comparisons=list(c("PFF", "PFF+LPS")), annotations="**",
              y_position = 0.27, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 0.34, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "PFF")), annotations="**",
              y_position = 0.41, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 0.5, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 0.5, tip_length = 0.02, size = 0.8, textsize=8, vjust=0.4) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,0.6)) +
  scale_fill_viridis_d()

plot.ens.pff.neuron.area

#Avg Int Neuron
mod.asyn.neuron.int <- lmer(asyn_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + 
                              (1|attempt.no) + (1|attempt.no:well), data = ens4)
warp.ens.asyn.avg.int.neuron <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no) + (1|attempt.no:well)
         , data = ens4))
summary(em.asyn.neuron.int <- emmeans(warp.ens.asyn.avg.int.neuron, ~Groups.f, type = "response"))
pwpm(em.asyn.neuron.int)

mod.pff.neuron.int <- lmer(pff_roi_tot_avg_int^(1/6) ~ 1 + Groups.f + (1|attempt.no), data = ens4)
warp.ens.pff.avg.int.neuron <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + 
        (1|attempt.no)
         , data = ens4))
summary(em.pff.neuron.int <- emmeans(warp.ens.pff.avg.int.neuron, ~Groups.f, type = "response"))
pwpm(em.pff.neuron.int)

plot.ens.asyn.neuron.int <- ggplot(ens4, aes(x=Groups.f, y=asyn_roi_tot_avg_int, fill=Groups.f)) +
  ylab("Average a-Syn Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Neuron a-Syn Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("6XPFF", "Control")), annotations="*",
              y_position = 93, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("LPS", "Control")), annotations="*",
              y_position = 85, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("LPS", "But+LPS")), annotations="**",
              y_position = 85, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="*",
              y_position = 85, tip_length = 0.05, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = asyn_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) + 
  coord_cartesian(ylim = c(0,100)) +
  scale_fill_viridis_d()

plot.ens.asyn.neuron.int

plot.ens.pff.neuron.int <- ggplot(ens4, aes(x=Groups.f, y=pff_roi_tot_avg_int, fill=Groups.f)) +
  ylab("Average PFF (ATTO-590) Intensity") + 
  xlab("Groups") + 
  labs(
    title = "Neuron PFF Average Intensity") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_signif(comparisons=list(c("PFF", "Control")), annotations="***",
              y_position = 55, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("Control", "2XPFF")), annotations="***",
              y_position = 70, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "4XPFF")), annotations="***",
              y_position = 105, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("2XPFF", "6XPFF")), annotations="***",
              y_position = 120, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("4XPFF", "4XPFF+But")), annotations="***",
              y_position = 105, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_signif(comparisons=list(c("PFF", "PFF+But+LPS")), annotations="***",
              y_position = 55, tip_length = 0.04, size = 0.8, textsize=8, vjust=0.4) +
  geom_dotplot(data = ens4, aes(x = Groups.f, y = pff_roi_tot_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",
        text = element_text(size = 14),
        axis.text.x = element_text(angle = 90),
        axis.title.x = element_blank()) +
  coord_cartesian(ylim = c(0,130)) +
  scale_fill_viridis_d()

plot.ens.pff.neuron.int

em.ens.asyn.neuron.int<- plot(em.asyn.neuron.int, 
                                type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Neuron a-Syn Average Intensity",
       x = "Estimated Marginal Mean (average intensity a-Syn)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())

em.ens.asyn.neuron.int

em.ens.pff.neuron.int<- plot(em.pff.neuron.int, type = "response", comparisons = TRUE, plotIT=FALSE) + 
  labs(title = "Neuron a-Syn Average Intensity",
       x = "Estimated Marginal Mean (Average Intensity PFF)",
       y = "Groups") +
  theme(text = element_text(size = 14),
          axis.title.y = element_blank())

em.ens.pff.neuron.int
```



```{r b3t transform}
#Transformation ENS
# hist(ens_tubb$avg.int)
# hist(ens_tubb$t.avg.int <- ens_tubb$avg.int^(1/6))
# skewness(ens_tubb$t.avg.int,na.rm = FALSE)
# 
# hist(ens_tubb$meant.avg.int <- abs(ens_tubb$avg.int - mean(ens_tubb$avg.int)))
# skewness(ens_tubb$t.avg.int,na.rm = FALSE)
# 
# bc.ens.tubb <- boxcox(ens_tubb$avg.int ~ ens_tubb$Groups.f)
# (lambda <- bc.ens.tubb$x[which.max(bc.ens.tubb$y)])
# ens_tubb$box.avg.int <- ((ens_tubb$avg.int^lambda-1)/lambda)
# 
# skewness(ens_tubb$box.avg.int,na.rm = FALSE)
# hist(ens_tubb$box.avg.int)
# 
# ens_tubb$beta.avg.int <- ens_tubb$avg.int / 255
# hist(ens_tubb$beta.avg.int)
# 
# #Cortex
# hist(cort_tubb$avg.int)
# hist(cort_tubb$t.avg.int <- cort_tubb$avg.int^(1/6))
# skewness(cort_tubb$t.avg.int,na.rm = FALSE)
```


```{r b3t analysis}

#ENS
# aov.21a <- kruskal.test(avg.int ~ Groups.f, data = ens_tubb)
# aov.21a
# np.d21 <- pairwise.wilcox.test(ens_tubb$avg.int, ens_tubb$Groups.f,p.adjust.method = "BH")
# np.d21

# 
# m.ens.tubb.1 <- lm(avg.int^(1/6) ~ Groups.f, 
#                      ens_tubb)
# m.ens.tubb.2 <- lmer(box.avg.int ~1 + Groups.f +
#                        (1|attempt.no),
#                      ens_tubb)
# m.ens.tubb.3 <- lmer(box.avg.int ~ Groups.f + 
#                        (1|well), 
#                      ens_tubb)
# m.ens.tubb.4 <- lmer(box.avg.int ~ Groups.f + 
#                        (1|attempt.no) + (1|attempt.no:well), 
#                      ens_tubb)
# lrtest(m.ens.tubb.1,m.ens.tubb.2,m.ens.tubb.3,m.ens.tubb.4)
# anova(m.ens.tubb.2,m.ens.tubb.3,m.ens.tubb.4)
# plot(m.ens.tubb.4)
# summary(em.ens.tub <- emmeans(m.ens.tubb.4, ~Groups.f, type = "response"))
# pwpm(em.ens.tub)
# 
# plot.ens.b3t.int.box <- ggplot(ens_tubb, aes(x=Groups.f, y=avg.int, fill=Groups.f)) + geom_boxplot()

plot.ens.b3t.int <- ggplot(ens_tubb, aes(x=Groups.f, y=tubb_avg_int, fill=Groups.f)) +
  ylab("Average B3T Intensity (unit8)") + 
  xlab("Groups") + 
  labs(
    title = "Average B3T Intensity (ENS)") +
  geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
    stat_summary(fun.data = mean_se,  
                 geom = "errorbar", position = 'dodge', colour = "black") + 
  geom_dotplot(data = ens_tubb, aes(x = Groups.f, y = tubb_avg_int),
               alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
               stackdir = "center", stackratio = 1) +
  theme(legend.position="none",axis.text.x = element_text(angle = 90))

plot.ens.b3t.int

ens_tubb.merge$bin.nuclei <- cut(ens_tubb.merge$nuclei.no, breaks = c(0,3,6,9,12,20,30,50))
ens_tubb.merge$bin.nuclei <- factor(ens_tubb.merge$bin.nuclei)

ens_tubb.merge$tubb_avg_int <- as.numeric(ens_tubb.merge$tubb_avg_int)
ens_tubb.merge$bin.tubb <- cut(ens_tubb.merge$tubb_avg_int, breaks = c(0,8,20))
ens_tubb.merge$bin.tubb <- factor(ens_tubb.merge$bin.tubb)

ens_tubb.simp <- subset(ens_tubb.merge, Groups.f == "Control")
ens_tubb.simp$tubb_avg_int <- as.numeric(ens_tubb.simp$tubb_avg_int)

plot.ens.tubb.asyn <- ggplot(ens_tubb.simp, aes(x=tubb_avg_int, y=asyn_roi_tot_avg_int, color = bin.nuclei)) +
  geom_point() +
  labs(
    title = "B3T Average Intensity vs. a-Syn intensity within whole neuron (Untreated Control)",
    color = "Nuclei Count") +
  xlab("B3T Average Intensity") +
  ylab("Average a-Syn Intensity") +
  scale_y_continuous(limits = c(0, 15), breaks = seq(0, 15, 3)) +
  scale_x_continuous(limits = c(0, 10), breaks = seq(0, 10, 2))

plot.ens.tubb.asyn

plot.ens.tubb.asyn2 <- ggplot(ens_tubb.simp, aes(x=tubb_avg_int, y=asyn_avg_int, color = bin.nuclei)) +
  geom_point() +
  labs(
    title = "B3T Average Intensity vs. a-Syn intensity within total image (Untreated Control)",
    color = "Nuclei Count") +
  xlab("B3T Average Intensity") +
  ylab("Average a-Syn Intensity") +
  scale_x_continuous(limits = c(0, 16), breaks = seq(0, 16, 4))

plot.ens.tubb.asyn2

ens_tubb.val <- subset(ens_tubb.merge, tubb_avg_int >= 0)
quantile(ens_tubb.val$asyn_avg_int, probs = 0.9)
ens_tubb.val$bin.asyn <- ifelse(ens_tubb.val$asyn_avg_int > 14.3454 ,1,0)

plot.ens.tubb.asyn3 <- ggplot(ens_tubb.merge, aes(x=tubb_avg_int, y=asyn_avg_int, color = bin.nuclei)) +
  geom_point() +
  labs(
    title = "B3T Average Intensity vs. a-Syn intensity within total image (All Samples)",
    color = "Nuclei Count") +
  xlab("B3T Average Intensity") +
  ylab("Average a-Syn Intensity") +
  scale_x_continuous(limits = c(0, 16), breaks = seq(0, 16, 4))

plot.ens.tubb.asyn4 <- ggplot(ens_tubb.val, aes(x=tubb_avg_int, y=asyn_avg_int, color = bin.nuclei)) +
  geom_point() +
  labs(
    title = "B3T Average Intensity vs. a-Syn intensity within total image (All Samples)",
    color = "Nuclei Count") +
  xlab("B3T Average Intensity") +
  ylab("Average a-Syn Intensity")

plot.ens.tubb.asyn5 <- ggplot(ens_tubb.val, aes(x=nuclei.no, y=asyn_avg_int, color = tubb_avg_int)) +
  geom_point() +
  labs(
    title = "Nuclei Count vs. a-Syn intensity within total image (All Samples)",
    color = "B3T intensity") +
  xlab("Nuclei Count") +
  ylab("Average a-Syn Intensity")

plot.ens.tubb.asyn5

m1 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f, data = ens_tubb.val)
m2 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int, data = ens_tubb.val)
m3 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + nuclei.no, data = ens_tubb.val)
m4 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int + nuclei.no, data = ens_tubb.val)

anova(m1,m2,m3,m4)
anova(m2,m4)

m2a <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int +
              (1|attempt), 
            data = ens_tubb.val)
m2b <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int +
              (1|well), 
            data = ens_tubb.val)
m2c <- lmer(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int +
              (1|attempt) + (1|attempt:well), 
            data = ens_tubb.val)

lrtest(m2,m2a,m2b,m2c)
anova(m2a,m2b,m2c)

#Analyze/Plot Whole Image
warp.m1 <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_avg_int) ~ Groups.f + tubb_avg_int +
        (1|attempt.no) + (1|attempt:well),
         data = ens_tubb.val))
summary(em.tub.asyn.raw.int <- emmeans(warp.m1, ~Groups.f, type = "response"))
pwpm(em.tub.asyn.raw.int)

plot(em.tub.asyn.raw.int, comparisons = TRUE, type = "response") +
  labs(title = "Whole Image a-Syn Average Intensity (B3T Adjusted)",
       y = "Groups")

warp.m1.pff <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_avg_int) ~ Groups.f + tubb_avg_int +
        (1|attempt.no) + (1|attempt:well),
         data = ens_tubb.val))
summary(em.tub.pff.raw.int <- emmeans(warp.m1.pff, ~Groups.f, type = "response"))
pwpm(em.tub.pff.raw.int)

plot(em.tub.pff.raw.int, comparisons = TRUE, type = "response") +
  labs(title = "Whole Image PFF Average Intensity (B3T Adjusted)",
       y = "Groups")

#Analyze/Plot Neuron
warp.m2 <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f + tubb_avg_int +
        (1|attempt.no) + (1|attempt:well),
         data = ens_tubb.merge))
summary(em.tub.asyn.neuron.int <- emmeans(warp.m2, ~Groups.f, type = "response"))
pwpm(em.tub.asyn.neuron.int)

plot(em.tub.asyn.neuron.int, comparisons = TRUE, type = "response") +
  labs(title = "Neuron a-Syn Average Intensity (B3T Adjusted)",
       y = "Groups")

warp.m2.pff <- with(fixed.avg.int.tran, 
    lmer(linkfun(pff_roi_tot_avg_int) ~ Groups.f + tubb_avg_int +
        (1|attempt.no) + (1|attempt:well),
         data = ens_tubb.merge))
summary(em.tub.pff.neuron.int <- emmeans(warp.m2.pff, ~Groups.f, type = "response"))
pwpm(em.tub.pff.neuron.int)

plot(em.tub.pff.neuron.int, comparisons = TRUE, type = "response") +
  labs(title = "Neuron PFF Average Intensity (B3T Adjusted)",
       y = "Groups")

ens_tubb.low <- subset(ens_tubb.val, tubb_avg_int < 18)

low.m1 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f, data = ens_tubb.low)
low.m2 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int, data = ens_tubb.low)
low.m3 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + bin.nuclei, data = ens_tubb.low)
low.m4 <- lm(asyn_roi_tot_avg_int^(1/6) ~ Groups.f + tubb_avg_int + bin.nuclei, data = ens_tubb.low)

warp.low.m1 <- with(fixed.avg.int.tran, 
    lmer(linkfun(asyn_roi_tot_avg_int) ~ Groups.f +
        (1|attempt.no)
         , data = ens_tubb.val))
summary(em.pff.neuron.int.low <- emmeans(warp.low.m1, ~Groups.f, type = "response"))
pwpm(em.pff.neuron.int.low)



# plot.em.ens.tubb <- plot(em.ens.tub, type = "response", comparisons = TRUE, plotIT=FALSE) +
#   labs(title = "Average B3T intensity (ENS)",
#        x = "Estimated Marginal Mean (average intensity B3T)",
#        y = "Groups")
# plot.em.ens.tubb
# 
# #Cortex
# m.cor.tubb.1 <- lm(sqrt(avg.int) ~ Groups.f, 
#                      cort_tubb)
# m.cor.tubb.2 <- lmer(sqrt(avg.int) ~1 + Groups.f +
#                        (1|attempt.no),
#                      cort_tubb)
# m.cor.tubb.3 <- lmer(sqrt(avg.int) ~ Groups.f + 
#                        (1|well), 
#                      cort_tubb)
# m.cor.tubb.4 <- lmer(sqrt(avg.int) ~ Groups.f + 
#                        (1|attempt.no) + (1|attempt.no:well), 
#                      cort_tubb)
# lrtest(m.cor.tubb.1,m.cor.tubb.2,m.cor.tubb.3,m.cor.tubb.4)
# anova(m.cor.tubb.2,m.cor.tubb.3,m.cor.tubb.4)
# plot(m.cor.tubb.4)
# summary(em.cor.tub <- emmeans(m.cor.tubb.4, ~Groups.f, type = "response"))
# pwpm(em.cor.tub)

# plot.cort.b3t.int <- ggplot(cort_tubb, aes(x=Groups.f, y=avg.int, fill=Groups.f)) +
#   ylab("Average B3T Intensity (unit8)") + 
#   xlab("Groups") + 
#   labs(
#     title = "Average B3T Intensity (Cortex)") +
#   geom_bar(position = 'dodge', stat = 'summary', fun = 'mean') +
#   geom_signif(comparisons=list(c("4XPFF", "Control")), annotations="***",
#               y_position = 25, tip_length = 0.05, vjust=0.4) +
#   geom_signif(comparisons=list(c("PFF", "4XPFF")), annotations="***",
#               y_position = 20, tip_length = 0.05, vjust=0.4) +
#     stat_summary(fun.data = mean_se,  
#                  geom = "errorbar", position = 'dodge', colour = "black") + 
#   geom_dotplot(data = cort_tubb, aes(x = Groups.f, y = avg.int),
#                alpha = .8, position_dodge(width = 0.9), binaxis = "y", dotsize = 0.1,
#                stackdir = "center", stackratio = 1) +
#   theme(legend.position="none")
# 
# plot.cort.b3t.int
# 
# plot.em.ens.tubb <- plot(em.cor.tub, type = "response", comparisons = TRUE, plotIT=FALSE) +
#   labs(title = "Average B3T intensity (Cortex)",
#        x = "Estimated Marginal Mean (average intensity B3T)",
#        y = "Groups")
# plot.em.ens.tubb
```



